// -*-c++-*-

#ifndef _NR_ARRAY_H_
#define _NR_ARRAY_H_

template <class T> class NRVec {
 private:
  int nn;	// size of array. upper index is nn-1
  T *v;
 public:
  NRVec();
  explicit NRVec(int n);		// Zero-based array
  NRVec(const T &a, int n);	//initialize to constant value
  NRVec(const T *a, int n);	// Initialize to array
  NRVec(const NRVec &rhs);	// Copy constructor
  NRVec & operator=(const NRVec &rhs);	//assignment
  NRVec & operator=(const T &a);	//assign a to every element
  inline T & operator[](const int i) { return v[i]; } 	//i'th element
  inline const T & operator[](const int i) const { return v[i]; }
  inline int size() const {  return nn; }
  inline const T* array() const { return v; }
  ~NRVec();
};

template <class T>
NRVec<T>::NRVec() : nn(0), v(0) {}

template <class T>
NRVec<T>::NRVec(int n) : nn(n), v(new T[n]) {}

template <class T>
NRVec<T>::NRVec(const T& a, int n) : nn(n), v(new T[n])
{
  for(int i=0; i<n; i++)
    v[i] = a;
}

template <class T>
NRVec<T>::NRVec(const T *a, int n) : nn(n), v(new T[n])
{
  for(int i=0; i<n; i++)
    v[i] = *a++;
}

template <class T>
NRVec<T>::NRVec(const NRVec<T> &rhs) : nn(rhs.nn), v(new T[nn])
{
  for(int i=0; i<nn; i++)
    v[i] = rhs[i];
}

template <class T>
NRVec<T> & NRVec<T>::operator=(const NRVec<T> &rhs)
  // postcondition: normal assignment via copying has been performed;
  //		if vector and rhs were different sizes, vector
  //		has been resized to match the size of rhs
{
  if (this != &rhs)
    {
      if (nn != rhs.nn) {
	if (v != 0) delete [] (v);
	nn=rhs.nn;
	v= new T[nn];
      }
      for (int i=0; i<nn; i++)
	v[i]=rhs[i];
    }
  return *this;
}

template <class T>
NRVec<T> & NRVec<T>::operator=(const T &a)	//assign a to every element
{
  for (int i=0; i<nn; i++)
    v[i]=a;
  return *this;
}

template <class T>
NRVec<T>::~NRVec()
{
  if (v != 0)
    delete[] (v);
}

template <class T>
class NRMat {
 private:
  int nn; // rows
  int mm; // columns
  T **v;
 public:
  NRMat();
  NRMat(int n, int m);			// Zero-based array
  NRMat(const T &a, int n, int m);	//Initialize to constant
  NRMat(const T *a, int n, int m);	// Initialize to array
  NRMat(const NRMat &rhs);		// Copy constructor
  NRMat & operator=(const NRMat &rhs);	//assignment
  NRMat & operator=(const T &a);		//assign a to every element
  inline T* operator[](const int i);	//subscripting: pointer to row i
  inline const T* operator[](const int i) const;
  inline int nrows() const;
  inline int ncols() const;
  ~NRMat();
};

template <class T>
NRMat<T>::NRMat() : nn(0), mm(0), v(0) {}

template <class T>
NRMat<T>::NRMat(int n, int m) : nn(n), mm(m), v(new T*[n])
{
  v[0] = new T[m*n];
  for (int i=1; i< n; i++)
    v[i] = v[i-1] + m;
}

template <class T>
NRMat<T>::NRMat(const T &a, int n, int m) : nn(n), mm(m), v(new T*[n])
{
  int i,j;
  v[0] = new T[m*n];
  for (i=1; i< n; i++)
    v[i] = v[i-1] + m;
  for (i=0; i< n; i++)
    for (j=0; j<m; j++)
      v[i][j] = a;
}

template <class T>
NRMat<T>::NRMat(const T *a, int n, int m) : nn(n), mm(m), v(new T*[n])
{
  int i,j;
  v[0] = new T[m*n];
  for (i=1; i< n; i++)
    v[i] = v[i-1] + m;
  for (i=0; i< n; i++)
    for (j=0; j<m; j++)
      v[i][j] = *a++;
}

template <class T>
NRMat<T>::NRMat(const NRMat &rhs) : nn(rhs.nn), mm(rhs.mm), v(new T*[nn])
{
  int i,j;
  v[0] = new T[mm*nn];
  for (i=1; i< nn; i++)
    v[i] = v[i-1] + mm;
  for (i=0; i< nn; i++)
    for (j=0; j<mm; j++)
      v[i][j] = rhs[i][j];
}

template <class T>
NRMat<T> & NRMat<T>::operator=(const NRMat<T> &rhs)
  // postcondition: normal assignment via copying has been performed;
  //		if matrix and rhs were different sizes, matrix
  //		has been resized to match the size of rhs
{
  if (this != &rhs) {
    int i,j;
    if (nn != rhs.nn || mm != rhs.mm) {
      if (v != 0) {
	delete[] (v[0]);
	delete[] (v);
      }
      nn=rhs.nn;
      mm=rhs.mm;
      v = new T*[nn];
      v[0] = new T[mm*nn];
    }
    for (i=1; i< nn; i++)
      v[i] = v[i-1] + mm;
    for (i=0; i< nn; i++)
      for (j=0; j<mm; j++)
	v[i][j] = rhs[i][j];
  }
  return *this;
}

template <class T>
NRMat<T> & NRMat<T>::operator=(const T &a)	//assign a to every element
{
  for (int i=0; i< nn; i++)
    for (int j=0; j<mm; j++)
      v[i][j] = a;
  return *this;
}

template <class T>
inline T* NRMat<T>::operator[](const int i)	//subscripting: pointer to row i
{
  return v[i];
}

template <class T>
inline const T* NRMat<T>::operator[](const int i) const
{
  return v[i];
}

template <class T>
inline int NRMat<T>::nrows() const
{
  return nn;
}

template <class T>
inline int NRMat<T>::ncols() const
{
  return mm;
}

template <class T>
NRMat<T>::~NRMat()
{
  if (v != 0) {
    delete[] (v[0]);
    delete[] (v);
  }
}

template <class T>
class NRMat3d {
 private:
  int nn;
  int mm;
  int kk;
  T ***v;
 public:
  NRMat3d();
  NRMat3d(int n, int m, int k);
  inline T** operator[](const int i);	//subscripting: pointer to row i
  inline const T* const * operator[](const int i) const;
  inline int dim1() const;
  inline int dim2() const;
  inline int dim3() const;
  ~NRMat3d();
};

template <class T>
NRMat3d<T>::NRMat3d(): nn(0), mm(0), kk(0), v(0) {}

template <class T>
NRMat3d<T>::NRMat3d(int n, int m, int k) : nn(n), mm(m), kk(k), v(new T**[n])
{
  int i,j;
  v[0] = new T*[n*m];
  v[0][0] = new T[n*m*k];
  for(j=1; j<m; j++)
    v[0][j] = v[0][j-1] + k;
  for(i=1; i<n; i++) {
    v[i] = v[i-1] + m;
    v[i][0] = v[i-1][0] + m*k;
    for(j=1; j<m; j++)
      v[i][j] = v[i][j-1] + k;
  }
}

template <class T>
inline T** NRMat3d<T>::operator[](const int i) //subscripting: pointer to row i
{
  return v[i];
}

template <class T>
inline const T* const * NRMat3d<T>::operator[](const int i) const
{
  return v[i];
}

template <class T>
inline int NRMat3d<T>::dim1() const
{
  return nn;
}

template <class T>
inline int NRMat3d<T>::dim2() const
{
  return mm;
}

template <class T>
inline int NRMat3d<T>::dim3() const
{
  return kk;
}

template <class T>
NRMat3d<T>::~NRMat3d()
{
  if (v != 0) {
    delete[] (v[0][0]);
    delete[] (v[0]);
    delete[] (v);
  }
}

#endif
